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1. Introduction 

QCD with three or more massless quark flavors would have a first order finite temperature 
transition [1], With increasing quark masses, the first order transition ends at a second order critical 
point, where the correlation length of the flavor singlet scalar meson diverges [2], 

Lattice QCD provided evidences [3, 4] of small singlet scalar meson screening masses, by 
simulations of the naive staggered fermion action, using a temporal extend, N t = 4. Since then, 
no new simulations have been exploring this area of research, partly because of the increasing 
difficulties of locating first order thermal transitions with improved staggered fermion formulations. 

Simulating 3-flavor finite temperature lattices using the nonperturbatively improved clover 
fermion action and the Iwasaki gauge action, we are able to locate the first order transition with 
lattice sizes of 24 3 x 8 and 32 3 x 8. This gives us the opportunity of studying the screening masses 
near the critical endpoint on the first order transition line. 

2. Measurement techniques 

We measure the connected screening propagators on three spatial directions, using multiple 
point sources either with or without Gaussian smearing, and point sinks. We use Gaussian smearing 
parameters, r = 4 and n = 100, so the source vector centered at xo is 

/ r 2 \" 

(j)(x)= n+—V 2 J S(x-x 0 ), (2.1) 

where 

V 2 *(x) = L Ui(x)x(x + i) + uj(x-t)x(x-i)-2x(x). (2.2) 

(= 1 , 2,3 

The flavor singlet screening correlator after zero momentum projection for Nf quark flavors, 

<£[«?(*)«?(*)] £[?Cy)</G0]> = Vf(W) 2 +A{e~ m + ■■■ (2.3) 

q q 

= -Nf(tr[S(x,y)rf(y,x)]) + N}(tr[S(x,x)]tr[S(y,y)]), (2.4) 

has a nonzero vacuum state due to finite chiral condensate, and requires evaluating disconnected 
diagrams. We call the exponential decay part the screening propagator, and write it as 

e -m\x\ +e -m\N l -x\ x _ c(jc) +NfD ( X ), (2.5) 

C(x) = (tr[5(x 0 ,xo-t-x)5 t (xo + x,x 0 )]), (2.6) 

D(x) = (tr[S(x 0 ,xo)]tr[S(xo+x,xo + x)]) - (tr[S(x 0 ,x 0 )]) 2 (2.7) 

= ^(tr[S(x 0 ,x 0 )] - {W)) 2 ) - (^( tr [S(xo,x 0 )] -tr[5(x 0 +x,x 0 +x)]) 2 ^) (2.8) 

where averaging over xo is implicit for brevity, Ni is the spatial extent of the lattice, S is the unsmeared 
quark propagator with ^tr[S(xo,xo)]> = and C and D are the connected and disconnected 

parts of the propagator respectively. The separation of the two terms in Eq. 2.8 makes them sensitive 
to the statistical noise from evaluating tr[S(x,x)]. They help us determine the required statistics and 
the computational cost for our method detailed below. 
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^min 

No. of diluted vectors 


12 3 4 

1 2 23 16 

1.14 0.82 0.58 


5 6 7 8 

120 210 411 256 

0.28 


Table 1: Number of space-time diluted vectors for a 32 3 x 8 lattice corresponding to different values 
of d m in, and the resulting improvement ratio of standard deviations. a imp is the standard deviation of 
(YW) computed with spin-color separation and space-time dilution; a is from spin-color combined 
undiluted stochastic sources. Standard deviations are estimated at the same computational cost at a 
fixed CG iteration number of 25. 


We use complex Z( 2)<r x Z(2)3 volume random sources, applying the idea of probing—instead 
of diluted uniform sources in the reference [5], we use stochastic sources with space-time dilution 
using the greedy multi-coloring algorithm. We test a few configurations with different values of 
r/ m irn which is the minimum number of links between nonzero sites in a diluted vector. We also 
combine it with spin-color separation [6], for the extra benefit verified in our tests. Table 1 shows 
a representative test result, for a particular gauge configuration of a lattice size 32 3 x 8. While 
increasing d m j n , the number of space-time diluted vectors increases exponentially with the exception 
of d m i n = 2, 4, and 8, which give relatively small number of diluted vectors due to symmetry. 
For d m i n = 1, which means no space-time dilution, spin-color separation alone does not help in 
reducing the standard deviation for the same cost. Using spin-color separation and space-time 
dilution together, we find increasingly smaller standard deviations with increasing d m \ n at the same 
computational cost. 

The improvement factor of four (with d min = 8 in table 1) comes if we solve 256 x 12 vectors 
for one configuration, which is impractical. The truncated solver method (TSM) [7] brings us further 
improvements when we extend it to use multiple levels of conjugate gradient (CG) truncations and 
apply spin-color separation and space-time dilution with different <7 mm at different truncation levels 
hierarchically. This is our method of hierarchical truncations with stochastic probing (HTwSP). 

One possible caveat: at a fixed CG iteration before converging, given the same initial guess, 
space-time dilutions with different d m \ n produce different solutions. Combining results from different 
ri m i n at different CG truncation levels would potentially contain bias. 

To test how big the difference is, we compare the undiluted results from 49152 volume random 
sources, and the diluted results from 256 source vectors per spin-color with space-time dilution of 
rimin = 4, which is also 49152 vectors (= 256 x 12 x 16). Figure 1 shows (\j/ t/r) at each spatial slices 
for one tested configuration of an ordered-start ensemble with a lattice size of 32 3 x 8, K = 0.14085, 
and /3 = 1.73. The figure shows results from the first eight CG iterations with a zero vector as the 
initial guess. Results from these two schemes are different at finite CG iterations. This difference, 
however, appears to decay exponentially. At a few hundred iterations where truncated results are 
used in this study, though untested, the difference should be well below statistical errors or even 
machine precision, so that no bias can be introduced in this way. 

For the HTwSP scheme in this study, we pick a handful of gauge configurations and search for 
the best combinations of CG truncation levels and values of d m { a for each level, while respecting our 
cost constraints and statistical precision needed for our study. Table 2 shows the best parameters 
we find on one of the test configurations, for the original one level TSM, two level and three level 
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(b) Difference as a function of CG iteration 


Figure 1: Results from the first eight CG iterations using undiluted random sources and diluted 
random sources of spin-color separation and space-time dilution with d m ; n = 4. 



(a) P = 1.73, K = 0.14085 (b) p = 1.74, k = 0.14054 

Figure 2: Lower panels show the connected and disconnected screening propagators of the flavor 
singlet scalar meson. Upper panels show the effective mass obtained from the full propagators, and 
the ones from disconnected propagators without error bars using light color as comparison. In the 
keys, ‘dis’ stands for the disordered-start ensemble, and ‘ord’ stands for the ordered-start ensemble. 


TSM (2L-TSM, 3L-TSM), and two level HTwSR We try to measure and minimize the computed 
ratios of variances, 7\' mip = G? /o 2 , improved against the undiluted one. For a fixed cost equivalent 
to about 500 convergent CG inversions, we find more than 40 times speedup using the scenario of 
2L-HTwSP. Further speedup is guaranteed with a larger cost constraint. 

3. Simulation results 

We start simulations with ordered and disordered gauge configurations with lattice sizes of 
24 3 x 8 and 32 3 x 8. For our two lattice coupling values of /3 = 1.73 and 1.74, we have tested a few 
K values and located the two-states signal at K = 0.14085 and 0.14054 respectively. Table 3 shows 
simple ensemble statistics. We measure screening propagators every 20 trajectories. 

Figure 2 shows the connected and disconnected parts of the flavor singlet scalar meson propa¬ 
gator and its effective masses as a solution to propagator values at nearby two spatial separations. 
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Scenario 

TSM 

2L-TSM 

3L-TSM 

2L-HTwSP 

Cost 

500.3 

493.2 

486.2 

501.9 

^imp 

0.203(4) 

0.119(4) 

0.102(3) 

0.0231(6) 

Ah 

99 

21 

5 

1 

Itei'h 

2049 

2049 

2049 

2049 

c h 

1 

1 

1 

24 

A 

2990 

1312 

442 

4 

ItCI'l 

275 

425 

475 

475 

Ci 

1 

1 

1 

192 

Al2 


8200 

2964 

1 

Iter 12 


50 

150 

200 

Cl2 


1 

1 

3072 


A B 13255 

Hero 25 

C 13 1 


Table 2: Improvements of different number of truncation levels, and one particular HTwSP scheme. 
Parameters are constrained by the computational cost. A ‘Cost’ of 1 equals computational cost of 
one convergent CG inversion. R is the ratio of variances (improved against undiluted one), which 
equals the inverse of the speedup in the wall-clock time for a fixed statistical error. While the 
subscripts denote different truncation levels, parameters for each level are: N is the number of 
random vector sets, the number of undiluted stochastic sources; C is the number of random vectors 
in one set, counting in dilution cost; Tter’ is the number of fixed CG iterations on that truncation 
level (2049 is the convergent iteration number). 


P 

K 

A/ 

start 

^conf 

plaquette 

gauge 

polyakov 

1.73 

0.14085 

24 

ord 

2930 

0.528741(67) 

1.24940(19) 

0.00831(16) 

1.73 

0.14085 

32 

ord 

13875 

0.528701(34) 

1.249514(98) 

0.008315(56) 

1.73 

0.14085 

24 

dis 

4060 

0.51495(15) 

1.28832(43) 

0.002555(73) 

1.73 

0.14085 

32 

dis 

16815 

0.514765(55) 

1.28883(16) 

0.002480(22) 

1.74 

0.14054 

32 

ord 

16000 

0.529832(44) 

1.24603(13) 

0.009558(32) 

1.74 

0.14054 

32 

dis 

10450 

0.520610(98) 

1.27195(28) 

0.003606(65) 


Table 3: The expectation values of plaquette, gauge density, and Polyakov loop. Number of 
configurations has first 1000 removed as thermalization. 


Because the disconnected part are orders of magnitude larger than the connected part, at long spatial 
separations, the effective masses from the full singlet propagator and the effective masses from the 
disconnected part alone agree. We see no statistical difference between the screening masses from 
disordered-start and ordered-start ensembles, which respectively correspond to the zero temperature 
and the finite temperature phases. 

Figure 3 shows meson screening masses at the transition. The flavor singlet scalar a meson 
is about half as light as the pion, and continues to decrease as the system approaches the critical 
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(a) j3 = 1.73, K = 0.14085 (b) p = 1.74, K = 0.14054 


Figure 3: Screening mass hierarchy across the first order transition. Symbols with light color 
represent results from 24 3 x 8 lattices, while the rest are from 32 3 x 8 lattices. 

endpoint, from /3 = 1.73 to j3 = 1.74. In the two parameter sets we have, the screening masses of 
the p, the n, and the a mesons are almost constant across the transition. We see the pion screening 
mass is smaller in the deconfined phase than the confined phase at /J = 1.73, K = 0.14085, while 
it is in an opposite order in the other parameter set. This is different from previous results [3, 4] 
with staggered fermions. As the system undergoes the first order transition from confined zero 
temperature phase to deconfined finite temperature phase, the a\ meson becomes degenerate with 
the p meson, and the ao meson drops and becomes closer to the n meson, as expected. 

Near the second order critical point, the system should exhibit the behavior of critical scalings. 
The flavor singlet scalar meson screening masses should scale as the inverse of the correlation 
length, m a °< t v , as proposed by Gavin et al. [2]. Other observables, we expect, should scale as 
a combination of the magnetic like order parameter, °< (we write it as fi to avoid confusion 
with the lattice coupling), and the energy like quantities, °< t~ a . Since we lack enough data sets 
to differentiate different universality classes, we treat the data as the mean field theory, such that 
a = 0 and j3 = v = 1/2. With the choice of the reduced temperature t = (j8 — jSc)//3c, where fic is 
the lattice coupling at the critical point, if the system is sufficiently close to the critical point, mfl 
should be linear in /J, and the squared discontinuities of the plaquette, the gauge density, and the 
chiral condensate (tshould also be linear in (5. Figure 4 shows these squared quantities, with 
lineal - lines going through two data points of each observable. Given the noisy scalar singlet meson 
screening masses, our extrapolations give consistent /3c- The mean field theory thus predicts that the 
critical endpoint is between /J = 1.745 and 1.75, which is consistent with our results from analysis 
of kurtosis intersections [8]. However, we need more data sets to determine the correct universality 
class. 

4. Summary 

We locate the first order transition with two-states signal at two parameter sets using the 
nonperturbatively improved clover fermion action and the Iwasaki gauge action. 

We developed an easy-to-implement method of hierarchical truncations with stochastic probing, 
which gives us an estimated 40 x speedup in measuring the quark condensate at spatial lattice slices. 

Our data show the a screening mass is about half of the pion mass on both sides of the 
first order transition close to the endpoint. A preliminary extrapolation to the critical endpoint 
shows consistency between our data and the mean field critical scaling. To determine the correct 
universality class requires more statistics and expanded parameter sets. 
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Figure 4: Linear extrapolations of m 2 a and the squared discontinuities of plaquette, gauge density, 
and (\j/y). 
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